Criticality of Vacancy-Induced Metal-Insulator Transition in Graphene 
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The criticality of vacancy-induced metal-insulator transition (MIT) in graphene is investigated by 
Kubo-Greenwood formula with tight-binding recursion method. The critical vacancy concentration 
for the MIT is determined to be x c « 0.053%. The scaling laws for transport properties near 
the critical point are examined showing several unconventional 2D localization behaviors. Our 
theoretical results have shed some new lights to the understanding of recent experiments in H-dosed 
graphene [Phys. Rev. Lett. 103, 056404 (2009)] and of 2D disordered systems in general. 
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Graphene, a single atomic layer of graphite, has drawn 
much interest in the last few years [l|. Because of the 
linear dispersion relation near the Dirac point, many un- 
conventional transport behaviors have been shown in the 
pristine 2D graphene (l-Q. On the other hand, defects 
and impurities are often present in graphene, so their 
effects on graphene transport properties are of signifi- 
cant scientific interests and practical implications. Va- 
cancy, one of the natural defects present in graphene, 
as observed in experiment yL has been a subject of 
intense theoretical study [6|4l2(. Vacancy is predicted 
to induce <masz- localized states in the vicinity of Fermi 
energy 0, H, G3 , leading to a metal-insulator transi- 
tion (MIT). However, some fundamental properties, es- 
pecially the criticality of the vacancy-induced MIT in 
graphene remains unclear. One important question is 
what is the critical vacancy concentration for MIT. Fur- 
thermore, what are the scaling laws characterizing the 
transport properties close to the critical point? The an- 
swers to these questions will help us not only to explain 
recent experiment results in graphene [HI but also to 
better understand, in general, localization and transport 
behaviors in 2D disordered systems. 

In this Letter, we report a theoretical study of crit- 
icality of the vacancy-induced MIT in graphene. Most 
importantly, we determine the critical vacancy concen- 
tration for MIT to be x c ~ 0.053%, which may be com- 
pared to the experimentally observed critical H coverage 
of ~ 0.08% for the H-dosing-induced MIT in graphene 
fl3j . We found that electrical conductivity a xx scales 
with the diffusion length L in a power-law, in contrast 
to the logarithmic scaling predicted for a conventional 
2D localization behavior 1 14, 

[H]. Above x r , there exists 
a mobility edge E c defining a transition point of energy 
from the localized states to the extended states, and the 
transition is discontinuous as characterized by a non-zero 
minimum conductivity o- m i n at E c . o- m in is found to vary 
from 1 to 1.3e 2 /ft, different from the universal value of 
0.628e 2 //z expected in a conventional 2D localization sys- 
tem [14J. Our results are consistent with the most salient 



features of the H-dosing-induced MIT in graphene ob- 
served recently [13], suggesting that both vacancy and H 
belong to short-range disorder in graphene. 

The graphene is described by the 7r-band, nearest- 
neighbor tight-binding Hamiltonian as 



(1) 



where i and j are the neighboring sites on the lattice, 
a — (t?i) is s P m index. S{ is the on-site energy set 
to zero, and 7^ = 7 = 2.7 eV is the hopping energy. 
To calculate the conductivity, the real-space recursion 
method [16] in Kubo-Greenwood formula is adopted, 
which has been proved to be a powerful tool in treat- 
ing various disordered graphene systems 17, [Hj]. The 
DC conductivity a xx at zero temperature can be calcu- 
lated as cr xx (E) = (e 2 /S)p(E) lim t _> +00 V(E, t). Here 
S is area per atom. p{E) = Tr[S(E — H)] is average 
density of states (ADOS). V(E,t) is diffusion coefficient, 

which is defined as V(E, t) = \ Tr{[£fi % S g:if ,&m , 
where x is the x-component of the position operator, 
U(t) = exp(— i^) is time evolution operator, and t is 
diffusion time. Periodic boundary conditions are used 
with a supercell sizes L x « 738 nm and L y « 213 nm. 
Vacancies are introduced simply by removing the atoms 
from the lattice 0,0]. 

Figure 1 shows the time dependence of diffusion coef- 
ficient V at Fermi energy (Ep = 0.0 eV) for different 
vacancy concentration (x) in the range of x > 0.1%. 
The inset shows the behavior of V at large time scale. 
Three main features are observed: (1) V first increases 
linearly with time and reaches its maximum V max at a 
very short time, less than 5 fs. V max is much smaller 
than V = Vpt for perfect graphene, where Fermi velocity 
v F = 3ta/2h ~ c/300. This means that the transport 
leaves the ballistic regime (dV/dt > 0) quickly and en- 
ters the diffusive regime (dV/dt = 0) due to vacancies. 
However, the diffusive regime is very short, without a 
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The relationship between conductivity and diffusion 
length is very useful to better understand the transport 
behavior. By definition, L(E,t) = y/T>(E, t)t, then we 

can easily find that L cx = t a , which is verified for all 
the concentrations simulated, as shown in Fig. 2(a) with 
a Log-Log plot of diffusion length (L) versus time (t). 

2a — 1 

Substituting t with L, <j xx scales with L as <j xx cx L^^~ . 
By dimensional analysis, (7 xx can be expressed as 



c (L)=a (e 2 /h)(L/iy 



(2) 



Figure 1: (Color online) Diffusion coefficient (V) as a function 
of time (£) for different vacancy concentration (x) at Ef- 
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Figure 2: (Color online) (a) Diffusion length (L) as a function 
of time (£) for different vacancy concentration (x) at Ef- (b) 
Scaling exponent a as a function of vacancy concentration (x) 
at Ef- 



plateau often seen in other systems. (2) Beyond V max , 
the system enters the localization regime, where V de- 
cays in a power law V cx t -A , similar to the behavior 
of disordered carbon nanotude [l9j]. This indicates that 
the system stays in the localization regime (dV/dt < 0) 
when x > 0.1% at large time scale, in agreement with 
the previous analytical results 0, 0, |l(j. (3) V decreases 
with the increasing x (inset of Fig. 1). This is because 
higher vacancy concentration leads to stronger scatter- 
ing, reducing the diffusion coefficient. We note that at a 
given vacancy concentration, p is a constant independent 
of time. Therefore, the conductivity a xx ~ Vp will have 
the same power-law scaling with time, a xx cx t -A , as V. 
Moreover, higher concentration will cause faster decay 
(larger A) in both diffusion coefficient and conductivity. 



where ao is the conductivity in the diffusive regime in 
unit of e 2 /h and I is the electron mean free path. It is im- 
portant to note that such power-law dependence of con- 
ductivity on diffusion length is different from the usual 
logarithmic dependence, <j(L) = ao — (e 2 /hn 2 Yln(L/l), as 
predicted by 2D scaling theory of localization [14]. We at- 
tribute this difference to the distinct nature of the quasi- 
localization induced by vacancies. The wave- function 
amplitude of the localized state induced by vacancies de- 
cays with the distance as 1/r [6, 8], which is not normal- 
izable in 2D. Therefore, such quasi-localized state [6, 8] 
is different from the usual localized state, whose wave 
function decays exponentially with the distance 

It has been predicted that the vacancy- induced quasi- 
localization will lead to MIT in graphene However, 
the critical vacancy concentration where the MIT occurs 
remains unknown. Here, we will determine the critical 
vacancy concentration for the MIT. From Fig. 2(a), we 
derive the scaling exponent, a, from the linear fitting of 
L versus t, and the resulting a are plotted as a function 
of x in Fig. 2(b). The curve of a(x) can be fit nicely by 
L cx t~^~ = t a , in the whole range of concentration. No- 
tice that when a = 0.5, A = 0, then dV/dt = because 
V cx t -A , which means a = 0.5 defines the position of dif- 
fusive regime, a > 0.5 (A < 0) means the ballistic regime 
and a < 0.5 (A > 0) means the localization regime. Thus, 
the power-law scaling holds true for all three transport 
regimes, a initially decreases rapidly with the increas- 
ing x, indicating the system quickly leaves the ballistic 
regime entering the localization regime. The transition 
point at a = 0.5 defines the critical vacancy concentra- 
tion for MIT, which is found to be x c w 0.053% in Fig. 
2(b). This value seems to agree quite well with the re- 
cent experimental result of critical concentration of MIT 
(0.08%) in H-dosed graphene [l3j]. It implies that similar 
transport behavior exists in these two disordered systems 
that both vacancy defect and H impurity create similar 
localization phenomena in graphene. In addition, the ex- 
periment suggested other unknown defects rather than H 
may induce MIT at even lower defect concentrations [l3[ , 
which might be attributed, at least partly, to vacancies 
according to our calculation. 

Next, we examine another criticality of the vacancy- 
induced MIT in graphene, the mobility edge (E c ), which 



3 




0.1 1 |@ 1 2 3 4 5 



X-X c (%) X (%) 



Figure 3: (Color online) (a), (b) and (c) Conductivity (a xx ) 
as a function of energy (E) for vacancy concentration x — 
0.5%, x = 1.0% and x = 5.0% at the time t\ when T> reaches 
its maximum (blue solid line) and £2 ~ 36.7 fs (orange dash 
line). "E" indicates the extended state and "L" indicates 
the localized state, (d) The mobility edge (E c ) as a function 
of x — x c . (e) The minimal conductivity (a m in) at E c as a 
function of x. 
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Figure 4: (Color online) (a) Mean free path (/) as a function of 
energy (E) for different vacancy concentration (x). (b) Mean 
free path (/) as a function of vacancy concentration (x) at Ef- 
The blue solid line is a fit to the data by / oc x -1 ^ 2 . 



defines the critical energy where the system changes from 
the localized state to the extended state above the criti- 
cal vacancy concentration x c . We only need to consider 
x > x ci because when x < x c the electronic states become 
extended over the whole energy range. Figures 3(a), (b) 
and (c) show the energy dependence of conductivity at 
two different times for three vacancy concentrations of 
0.5%, 1.0%, and 5.0%, respectively. One time is chosen 
at where V reaches its maximum, t\(D = Dm ax ), and the 
other is chosen at a much later time (£2 ~ 36.7 fs >> t\). 
At the low energies closed to Ef, o~ xx (ti) > cr xx (t2), but 
at the high energies away from Ep, cr xx (ti) < v xx {t2)- 
This is because the low-energy electronic states are lo- 
calized whose conductivity decreases with time; while the 
high-energy electronic states are extended whose conduc- 
tivity increases with time. This is in qualitative agree- 
ment with the previous analysis of the inverse participa- 
tion ratio (IPR) of electron wavefunction p, [8|. Then, 
the mobility edge E c can be obtained as the energy for 
°xx{ti) = c X x(t 2 ). In Figs. 3(a), (b) and (c), we see 
that the mobility edge E c increases with the increasing 
vacancy concentration, as the increasing vacancy concen- 
tration expands the energy range of localization when 
x > x c . 

In Fig. 3(d), we plot E c as a function of x, in a Log- 
Log scale, which clearly shows a power-law dependence 
of E c on x. Fitting the simulation data, we obtain the 
scaling relation of E c oc (x — XcY, with the scaling expo- 
nent P ~ 0.81. Importantly, a non-zero minimum con- 
ductivity o-min appears at the mobility edge. Figure 3(e) 
shows that o- min first decreases and then increases with 



x, but the reason for such trend is unclear. The existence 
of o-min is possibly due to the quasi-localized nature of 
electronic states induced by vacancy in graphene. It in- 
dicates that the transition from the localized states to 
the extended states at the mobility edge is discontinu- 
ous as Mott are; ued [13, suggesting the break-down of 
one-parameter scaling theory. Also, our simulated o- m in 
varies from 1 to 1.3e 2 //i, different from the universal 
value of 0.628e 2 //i as expected in conventional 2D lo- 
calization system [14]. 

Last, we examine the electron mean free path /, an 
important length scale characterizing the degree of elec- 
tron scattering by disorder. In Fig. 4(a), the en- 
ergy dependence of I is shown, which is obtained by 
1(E) = V max (E)/v(E) [lj, 0]. v(E) is the average 
wave-packet velocity at energy E, extracted from the dif- 
fusion coefficient in the ballistic regime at short times as 
V(E,t) « v(E) 2 t. It is interesting to observe a minimal 
mean free path at Ef- This distinguishes vacancies from 
other disorders in graphene that give rise to a maximal 
mean free path at E F [IzLElI]. The minimal mean free 
path means that the lifetime of electrons is the shortest 
near Ef in the localization regime, violating graphene's 
ordinary Fermi liquid behavior [21]. Such an unusual be- 
havior has also been experimentally observed in graphene 
dosed with a high concentration of H atoms [13]. 

Near Ef, the mean free path I is much less than the 
average distance between vacancies (d). This is because 
the quasi-localized states substantially reduce I in the 
low-energy region. For example, when x = 0.1%, I is 
found to be 1.6 nm, while d is in the order of 5 nm. In 
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Fig. 4(b), we plot I as a function of x, which shows that 
/ scales with x in a power law I = Ax 5 (5 = —0.5) at E F . 
Fitting the simulation data, we obtained the parameter 
A — 0.508 nm. This power-law scaling has been derived 
before from the full self-consistent Born approximation 
(FSB A) [7], but the FSB A is limited to extremely low 
concentration and can not describe the localization be- 
havior near Ep and give the quantitative values of I. 

Our quantitative results indicate that vacancies are 
short-range disorders similar to H atoms in graphene as 
suggested by a recent experiment [13] . The similar re- 
sults between vacancy defect and H impurity can be un- 
derstood by the dramatic reconstruction of the energy 
spectrum near the Dirac point in both cases, while the 
Fermi energy is preserved. As short-range local scatters, 
both vacancy and H lead to strong localization by sig- 
nificant scattering. We note that at low defect /impurity 
concentration, the mean free path I is shorter for vacancy 
(1.6 nm at x — 0.1%) than for H atom (~ 4 — 5 nm at 
x = 0.1%). This is because vacancy produces a much 
stronger local scattering center (creating an infinite bar- 
rier) than H (inducing a structural transition from sp 2 
to sp 3 hybridization). However, we find that the rate of 
I decreases faster with x for H (5 < — 1) [13| than that 
for vacancy (S — —0.5). This is possibly because the 
scattering center produced by vacancy is more localized 
limited to single atomic site, while the scattering center 
produced by H atom is relatively more extended as the 
sp 3 hybridization induces reconstruction in the neighbor- 
ing atoms around it. Consequently, at high concentra- 
tion, multiple scattering by H atoms becomes more and 
more important, leading to a faster decay of L A differ- 
ent case was shown in graphene substituted with boron 
or nitrogen, where the long-range scattering is induced 
by chemical disorder with 5 = — 1 [18]. 

In conclusion, the criticality of vacancy-induced MIT 
in graphene has been investigated with quantitative sim- 
ulations. The MIT is shown to occur at the critical va- 
cancy concentration x c « 0.053%, caused by the vacancy- 
induced gimsi-localization. Several unconventional criti- 
cal behaviors have been revealed in the defected graphene 
with vacancies, including the power-law rather than log- 
arithmic scaling relation between conductivity and diffu- 
sion length and the existence of a minimum conductiv- 
ity at the mobility edge that differs from the universal 
value of 0.628e 2 //z expected in conventional 2D localiza- 



tion systems. Our results for vacancies are consistent 
with the most salient features of the H-dosing-induced 
MIT in graphene as observed in recent experiments jl3j |. 
suggesting that both vacancy and H belong to the same 
class of short-range disorders in graphene with H dos- 
ing involving stronger multiple scattering. Our studies 
also shed some new lights to the general understanding 
of critical behaviors of 2D disordered systems. 
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